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We consider a parametrically driven Klein-Gordon system describing micro- and nano-devices, 
with integrated electrical and mechanical functionality. Using a multiscale expansion method we 
reduce the system to a discrete nonlinear Schrodinger equation. Analytical and numerical calcula- 
tions are performed to determine the existence and stability of fundamental bright and dark discrete 
solitons admitted by the Klein-Gordon system through the discrete Schrodinger equation. We show 
that a parametric driving can not only destabilize onsite bright solitons, but also stabilize intersite 
bright discrete solitons and onsite and intersite dark solitons. Most importantly, we show that there 
is a range of values of the driving coefficient for which dark solitons are stable, for any value of 
the coupUng constant, i.e. oscillatory instabilities are totally suppressed. Stability windows of all 
the fundamental solitons are presented and approximations to the onset of instability are derived 
using perturbation theory, with accompanying numerical results. Numerical integrations of the 
Klein-Gordon equation are performed, confirming the relevance of our analysis. 

PACS numbers: 05.45.-a, 63.20.Pw, 85.85.+j 



I. INTRODUCTION 

Current advances in the fabrication and control of elec- 
tromechanical systems on a micro and nanoscale bring 
many technological promises pj. These include efficient 
and highly sensitive sensors to detect stresses, vibrations 
and forces at the atomic level, to detect chemical signals, 
and to perform signal processing [2]. As a particular 
example, a nanoelectromechanical system (NEMS) can 
detect the mass of a single atom, due to its own very 
small mass [3, 4]. 

On a fundamental level, NEMS with high frequency 
will allow research on quantum mechanical effects. This 
is because NEMS, as a miniaturization of microelec- 
tromechanical systems (MEMS), can contain a macro- 
scopic number of atoms, yet still require quantum me- 
chanics for their proper description. Thus, NEMS can 
be considered as a natural playground for a study of 
mechanical systems at the quantum limit and quantum- 
to-classical transitions (see, e.g., Ref. 5 and references 
therein). 

Typically, nanoelectromechanical devices comprise an 
electronic device coupled to an extremely high frequency 
nanoresonator. A large number of arrays of MEMS and 
NEMS resonators have recently been fabricated experi- 
mentally (see, e.g., Ref. [6]). One direction of research 
in the study of such arrays has focused on intrinsic local- 
ized modes (ILMs) or discrete breathers. ILMs can be 
present due to parametric instabilities in an array of os- 
cillators [7] . ILMs in driven a rray s of MEMS have been 
observed experimentally @, 0, [l^- 

Motivated by a recent experiment of Buks and 
Roukes [il that succeeded in fabricating and exciting an 
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array of MEMS and measuring oscillations of the res- 
onators, here we consider the equation [TT] 



DA2ipn - [1 - H cos{2u;pt)]ipn ± (1) 



which governs the oscillation amplitude of such an array. 
Equation ([1]) is a simplified model of that discussed in 
Ref. [H, subject to an assumption that the piezoelectric 
parametric drive is applied directly to each oscillator [13] • 
The variable cpn represents the oscillation amplitude of 
the nth oscillator from its equilibrium position, D is a 
dc electrostatic nearest-neighbor coupling term, is a 
small ac component with frequency 2u;p responsible for 
the parametric driving, A2(pn = ^n+i "^^n^^n-i is the 
discrete Laplacian, the dot denotes the derivative with re- 
spect to t, and the 'plus' and 'minus' signs of the cubic 
term correspond to a 'softening' and 'stiffening' nonlin- 
ear ity, respectively. Here, we assume ideal oscillators, so 
there is no damping present. The creation, stability, and 
interactions of ILMs in ([1]) with low damping and in the 
strong-coupling limit, have been investigated in Ref. [ll|. 
Here, we extend that study to the case of small coupling 
parameter D. 

In performing our analysis of the governing equation 
we introduce a small parameter e <C 1, and assume 
that the following scalings hold: 



D 



We then expand each cpn in powers of e, with the leading- 
order term being of the form 

- e (V^n(T2, T3, . . . )e-^^° + V^JT2, T3, . . . )e^^°) , (2) 



where = e^T. Then the terms at 0{e'e 
yield the following equation for ipn (see Refs 
a related reduction method): 



m O 
Hfor 



±i3AT2 



±3|V^n|Vn, (3) 



where the dot now denotes the derivative with re- 
spect to T2. Correction terms in Eq. ([2]) are of order 



2 



0(^g±i(/c+i)To^^3g±i(/c-i)To)^ k A justification of 

this rotating wave type approximation can be obtained 
in, e.g., Ref. fl6l . 

Writing ?An(T2) = (t)n{T2)e^'^^^'^^\ we find that 
Eq. ([3]) becomes 

2 • _ ^ 

- = C'A20n T Hn T 7^n ± l^nT^n- (4) 

Then, taking T2 = |T, we find that the equation 
above becomes the parametric driven discrete nonhnear 
Schrodinger (DNLS) equation 

Z^n = -CA20n ± A0n ± 7</>n T |0n|^(/>n; (5) 

here the dot denotes the derivative with respect to T. 
The softening and stiffening nonhnearity of ([T]) corre- 
spond, respectively, to the so-called focusing and defo- 
cusing nonhnearity in the DNLS ([5]). 

In the absence of parametric driving, i.e., for 7 = 0, 
Eq. ([5]) is known to admit bright and dark solitons in the 
system with focusing and defocusing nonhnearity, respec- 
tively. Discrete bright solitons in such a system have been 
discussed before, e.g. in Refs. [13, [3, where it was 
shown that one-excited-site (onsite) solitons are stable 
and two-excited-site (inter site) solitons are unstable, for 
any coupling constant C. Undriven discrete dark solitons 
have also been examined [20, 21,^22', 23,^24]; it is known 
that intersite dark solitons are always unstable, for any 
C, and onsite solitons are stable only in a small window 
in C. Furthermore, an onsite dark soliton is unstable due 
to the presence of a quartet of complex eigenvalues, i.e., 
it suffers oscillatory instability. 

The parametrically driven DNLS ([5]) with a focusing 
nonhnearity and finite C has been considered briefiy in 
Ref. where it was shown that an onsite bright dis- 
crete soliton can be destabilized by parametric driving. 
Localized excitations of the continuous limit of the para- 
metrically driven DNLS, i.e. ([5]) with C ^ 00, have been 
considered by Barashenkov et al. in a different context 
of applications [2I, [13, [H, [H, [sO, [m, |32|. The same 
equation also applies to the study of Bose-Einstein con- 
densates, describing the so-called long bosonic Josephson 
junctions [33l, 1341] . 

In this paper, we consider ([1]) with either softening 
or stiffening nonlinear it ies, which admit bright or dark 
discrete solitons, respectively. The existence and stabil- 
ity of the fundamental onsite and intersite excitations 
are discussed through the reduced equation ([5]). Eq. ([5]) 
and a corresponding eigenvalue problem are solved nu- 
merically for a range of values of the coupling and driv- 
ing constants, C and 7, giving stability windows in the 
(C, 7) plane. Analytical approximations to the bound- 
aries of the numerically obtained stability windows are 
determined through a perturbation analysis for small C. 
From this analysis, we show, complementing the result 
of Ref. [25I, that parametric driving can stabilize intersite 
discrete bright solitons. We also show that parametric 
driving can even stabilize dark solitons, for any coupling 



constant C. These findings, which are obtained from the 
reduced equation ([5]), are then confirmed by direct nu- 
merical integrations of the original governing equation 

The present paper is organized as follows. In Sec. II we 
present the existence and stability analysis of onsite and 
intersite bright solitons. Analysis of dark solitons is pre- 
sented in Sec. III. Confirmation of this analysis, through 
numerical simulations of the Klein-Gordon system ([1]) , is 
given in Sec. IV. Finally, we give conclusions in Sec. V. 



II. BRIGHT SOLITONS IN THE FOCUSING 
DNLS 

In this section we first consider the existence and sta- 
bility of bright solitons in the focusing DNLS equation. 
For a static solution of ([5]) of the form (j)^ = where 
Un is real- valued and time-independent, it follows that 



CA2Un - + + 7^n = 0. 



(6) 



Once such discrete solitary-wave solutions of ([5]) have 
been found, their linear stability is determined by solv- 
ing a corresponding eigenvalue problem. To do so, we 
introduce the linearization ansatz 

where (5 <C 1, and substitute this into ([5]), to yield the 
following linearized equation at 0{5): 



-CA2en - 2\Un\ en - U^Cn + ACn + 7^X1- (7) 



Writing en{t) = r]n -\- i^n^ we then obtain from Eq. ([7j) 
the eigenvalue problem 



in 



where 



n 



H 



£+(C) 
-C-{C) 



(8) 



and the operators C-{C) and Cj^iC) are defined by 

C-{C) = -CA2 - (3^2 - A - 7), 
£+(C) ^ -Ci^2-{ul-K + i). 

The stability of the solution Un is then determined by the 
eigenvalues of ([8]). If we denote these eigenvalues by iuo^ 
then the solution Un is stable only when Im(a;) = for 
all eigenvalues uo. 

We note that, because ([8]) is linear, we may eliminate 
one of the eigenvectors, for instance to obtain an alter- 
native expression of the eigenvalue problem in the form 



(9) 



In view of the relation 1^ = cj^, it follows that a soliton 
is unstable if it has an eigenvalue with either ^7 < or 
lm{n) ^ 0. 
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A. Analytical calculations 

Analytical calculations of the existence and stability 
of discrete solitons can be carried out for small coupling 
constant C, using a perturbation analysis. This analysis 
exploits the exact solutions of (|6]) in the uncoupled limit 

(7 = 0, which we denote by Un = Un \ in which each Un^ 
must take one of the three values given by 



0, ±VA + 7. 



(10) 



Solutions of (|6|) for small C can then be calculated ana- 
lytically by writing 

In studying the stability problem, it is natural to also 
expand the eigenvector r]n and the eigenvalue ft in powers 
of C, as 

Upon substituting this expansion into Eq. (|9]) and iden- 
tifying coefficients of successive powers of the small pa- 
rameter C, we obtain from the equations at 0{1) and 
0{C) the results 



£+(0)/:_(0)-l](^) 



(0) = 



(0) 
n 5 



(11) 

(12) 



where 



/„ = -(A2+2u(°)u(i))£_(0)-£+(0)(A2+6wf 

(13) 

In the uncoupled limit, C = 0, the eigenvalue problem 
is thus simphfied to 



=£+(0)£_(0), 



(14) 



from which we conclude that there are two possible eigen- 
values, given by 



n 



(0) 



7^ ^ 



(0) 
E 



4(A + 7)7, 



(0) 







which correspond, respectively, to the solutions u 
(for all n) and Un^ = ±v^A + 7 (for all n). 

We begin by considering bright soliton solutions, for 
which ^ as n ^ ±00. This then implies that (for 



C 



0) the eigenvalue has infinite multiplicity; it 
generates a corresponding continuous spectrum (phonon 
band) for finite positive C. To investigate the significance 
of this continuous spectrum, we introduce a plane wave 
expansion 

from which one obtains the dispersion relation 

n = (2C(cos K.-1)-Af- 7^ (15) 



y<yth 



y>yth 
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FIG. 1: A sketch of the dynamics of the eigenvalues and the 
continuous spectrum of a stable onsite bright soliton in the 
(Re(Q),Im(Q)) plane. The arrows indicate the direction of 
movement as the coupling constant C increases. Note that 
a soliton is unstable if there is some Q with either Q < or 
lm{Q) / 0. 



This in turn shows that the continuous band lies between 



= — 7^, when hz = 0^ 



(16) 



and 



A^ - 7^ + 8C(A + 2C), when = tt. (17) 



From the continuous spectrum analysis above, it can 
be concluded that an instability can only be caused by 
the dynamics of discrete spectrum. 



1. Onsite bright solitons 

The existence and stability of a single excited state, i.e. 
an onsite bright soliton, in the presence of a parametric 
driving has been considered in Ref. ,25|. For small C, the 
soliton is given by [25] 

f VA T^ + C/y/AT^ + (9(C^), n = 0, 
Un = { C/VAT^ + 0{C^), n = -1, 1, 

[ O(C^), otherwise, 

(18) 

and its eigenvalue by 

= 4(A + 7)7 + 87C + 0{C^). (19) 

It was shown in Ref. [25] that the configuration (p!8|) . 
which is known to be stable for any value of C when 
7 = 0, can be destabilized by parametric driving. Fur- 
thermore, it was shown that there are two mechanisms of 
destabilization, as sketched in Fig. [TJ The two instabil- 
ity scenarios are determined by the relative positions of 



r^^^ and ^^^^ as we now summarize. First, we note that 
there is a threshold value, 7th = A/5, at which the two 
leading-order eigenvalues coincide, so that 1]^^ = ^c^- 



4 



For 7 > 7th, upon increasing C from C = 0, the in- 
stability is caused by the cohision of VLe with Vtu] tak- 
ing VLe = then yields the corresponding approximate 
critical value 

7cr = -^A - + i ^ 56^A + 96C^. (20) 



For 7 < 7th, by contrast, the instability is caused by the 
collision of Qe with an eigenvalue bifurcating from Ql. 
In this case, the critical value of 7 can be approximated 
by taking Qe = ^L, giving 

7c'r = -^A-^C+ JV9A' + 16C(A + C). (21) 

Together, 7^^. and 7^^ give approximate boundaries of the 
instability region in the (C, 7)-plane. 

y<yth 



£2^ 



y>yth 




FIG. 2: As Fig. [T] but for a stable intersite bright soliton. 



2. Intersite bright solitons 

The next natural fundamental solution to be consid- 
ered is an intersite bright soliton, i.e., a two-excited- 
site discrete mode. In the uncoupled limit, the mode 
structure is of the form Un^ = for n 7^ 0, 1 and 
u"^^ = = \/A + 7. Using a perturbative expansion, 
one can show further that the soliton is given by 

f y/A T^ + iC/^AT^ + 0(C^), n = 0,l, 
Un = { C/VAT^ + (9(C2), n = -1, 2, 

[ 0{C^), otherwise. 

(22) 

To study the stability of the intersite bright soliton 
above, let us consider the 0{1) equation (pTj) . Due to 
the presence of two non-zero excited sites at C = 0, the 
soliton ([22]) has at leading order the double eigenvalue 

1^^^ = 4(A + 7)7, with corresponding eigenvectors r]n-' = 

for n 7^ 0, 1, r]Q^^ = 1, and r]^^^ = ±1. 

The continuation of the eigenvalue 1]^^ above for 
nonzero coupling C can be obtained from Eq. (p!2|) by ap- 
plying a solvability condition. The Fredholm alternative 



requires that /n = for all n, from which we immediately 
deduce that the double eigenvalue splits into two distinct 
eigenvalues, which are given as functions of C by 

= 4(A + 7)7 + 47C + (9(C'), (23) 

and 

^2 = 4(A + 7)7 - 4(A + 7)C + 0{C^). (24) 

As is the case for onsite discrete solitons, intersite 
bright solitons can also become unstable. The mecha- 
nism of the instability is again determined by the relative 
positions of 1]^^ and Q^^-' , as sketched in Fig. O Per- 
forming an analysis corresponding to that in Ref. [25[, 
we find that the two mechanisms of destabilization for 
an onsite discrete soliton also occur here. The two sce- 
narios have corresponding critical values of 7, which are 
given as functions of C by 

7cr = -7A + + iV9A2 + 52AC + 84C2, (25) 
5 5 5 

7cr = -^A - + iV9A2 + 8AC + 4C2. (26) 
000 

We emphasize, as is apparent from the sketch shown 
in Fig. [21 that there is another possible mechanism of 
destabilization for 7 < 7th, namely when becomes 
negative. The third critical choice of parameter values is 
then obtained by setting Q2 = 0, i.e. 

li = C. (27) 

B. Comparisons with numerical calculations 

We have solved the steady-state equation (|6]) numeri- 
cally using a Newton-Raphson method, and analyzed the 
stability of the numerical solution by solving the eigen- 
value problem ([8]). In this section, we compare these 
numerical results with the analytical calculations of the 
previous section. For the sake of simplicity, we set A = 1 
in all the illustrative examples. 

1. Onsite bright solitons 

Comparisons between numerical calculations and ana- 
lytical approximations for the case of onsite bright soli- 
tons have been fully presented and discussed in Ref. [25[ . 
For the sake of completeness, we reproduce the results of 
Ref. [25] for the (in) stability domain of onsite solitons in 
the (C, 7) plane in Fig. [3l Approximations ([2Q|) and ([2T]) 
are also shown there. 



2. Intersite bright solitons 

For the stability of intersite bright solitons, we start 
by examining the validity of our analytical prediction for 
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FIG. 3: (Colour online) The (in)stability region of onsite 
bright solitons in (C, 7) space. For each value of C and 7, the 
corresponding colour indicates the maximum value of |Im(cj)| 
(over all eigenvalues cj) for the steady-state solution at that 
point. Stability is therefore indicated by the region in which 
Im(a;) = 0, namely the dark region. White dashed and dash- 
dotted lines give the analytical approximations (|2Q|) and ()2ip , 
respectively. 
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FIG. 4: Comparison between the critical eigenvalue of inter- 
site bright solitons obtained numerically (solid lines) and its 
analytical approximation (dashed lines). The upper and lower 
curves correspond, respectively, to 7 = 0.5 and 7 = 0.1, ap- 
proximated by Eq. (|24)) . whereas the middle one corresponds 
to 7 = 0.18, approximated by Eq. (|23)) . 



the eigenvalue associated with the phase mode as given 
by Eqs. ([23|) and ([24|) . In Fig. [H we present a compar- 
ison between the analytical approximation and the nu- 
merics for some representative values of 7 (specifically 
7 = 0.1,0.18,0.5). This figure reveals the relative ac- 
curacy of the small-C approximations, and we conclude 
that their range of validity is wider for smaller values 
of 7. 

Next we turn to a description of the eigenvalue struc- 
ture of this intersite configuration for the three values of 
7 given above; this is shown in Fig. [5l where the left and 
right panels respectively present the structure just before 



and just after the first collision that results in the mode 
instability. We now describe results in more detail for 
the three values of 7 in turn. 

For 7 = 0.1, when C = the eigenvalues uj lie in the 
gap between the two parts of the continuous spectrum, 
and the instability is caused by a collision between the 
critical eigenvalue and its twin at the origin (see the top 
panels of Fig. [5|). For 7 = 0.18, the eigenvalues uo also 
lie in the gap between the two parts of the continuous 
spectrum, but the instability in this case is due to a col- 
lision between one of the eigenvalues and the inner edge 
of the continuous spectrum at a; = iV^H^ (see the mid- 
dle panels of Fig. [5]). In contrast to the two cases above, 
for 7 = 0.5 the eigenvalues lie beyond the continuous 
spectrum, and the instability is caused by a collision be- 
tween the critical eigenvalue and the outer boundary at 
uo = ±\^ftu (see the bottom panels of Fig. [5]). All the 
numerical results presented here are in accordance with 
the sketch shown in Fig. [2l 

Numerical calculations of the stability of intersite 
bright solitons, for a relatively large range of C and 7, 
give us the stability domain of the bright solitons in the 
two-parameter (C, 7) plane, which is presented in Fig. [6l 
We use colours to represent the maximum of |Im(a;)| as a 
function of C and 7; thus solitons are stable in the black 
region. Our analytical predictions for the occurrence of 
instability, given by Eqs. ([25]) -([27 |) . are also shown, re- 
spectively, by dashed, dotted, and dash-dotted lines. 



III. DARK SOLITONS IN THE DEFOCUSING 
DNLS 

In this section we consider the existence and stability of 
onsite and intersite dark solitons for the defocusing DNLS 
equation. Then a static (real-valued, time-independent) 
solution Un of ([5]) satisfies 



AUn - ^Un 0. 



(28) 
±00, 



In contrast to bright solitons, where ^ as n 
dark solitons have Un ±a/A 7 as n ±00. 

To examine the stability of Un^ we again introduce the 
linearization ansatz <pri = Un -\- Scn^ where again 5 <^ 1. 
Substituting this ansatz into the defocusing equation ([5]), 
writing en{t) = rjn -\- and linearizing in ^, we again 
find 



Vn 
in 











Vn 
in 



= n 



Vn 
in 



(29) 



but where the operators C±{C) are now defined as 

C-{C) = -CA2 + (3?|2 - A - 7), 
£+(C) = -CA2 + (u2 -A + 7). 

The eigenvalue problem above can be simplified further 
as for the focusing case, to the alternative form 



(30) 
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(a)7 = 0.1, C = 0.05 



(b)7 = 0.1, C = 0.3 
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Re((o) 

(c)7 = 0.18, C = 0.05 




-0.5 0.5 

Re(o)) 

(d)7 = 0.18, C = 0.18 





(e)7 = 0.5, C = 0.05 



(f)7 = 0.5, C = 0.2 



FIG. 5: The structure of the eigenvalues of intersite bright sohtons in the complex plane for three values of 7, as indicated in 
the caption of each plot. Left and right panels depict the eigenvalues of stable and unstable solitons, respectively. 



Performing a stability analysis as before, we find the dis- and so the continuous band lies between 

persion relation for a dark soliton to be ^ w * x ^ 

CIl = 4(A + 7)7, when k = 0, 

and 



(32) 



n = (2C{cos K-1)-(A + 2j) f - A\ (31) = 4(A + 7)7 + 8C(A + 27 + 2C), when k = tt. (33) 



7 



^- 0.2 




FIG. 6: (Colour online) As Fig. O but for intersite bright 
solitons. Our analytical approximations, given by Eqs. ()25p . 
([26|) and ([27|) . are shown as white dash- dotted, dotted, and 
dashed lines, respectively. 



A. Analytical calculations 

To study the eigenvalue (s) of the dark soliton analyt- 
ically, we again expand r]n and ft in powers of C, and 
hence obtain from (|30]), at 0{C^) and 0{C^), respec- 
tively, the equations 



and 



with 



where 



/:+(o)/:_(o)-f2(°)l7?f=o, 
i:+(o)£_(o) 



fn 



(34) 

(35) 
(36) 



Qn = (A2 - 2«(°)«(i))£_(0) + £+(0)(A2 - 6«(°)^^(l)). 

(37) 

We next investigate the eigenvalues of both intersite and 
onsite modes. 



From Eq. ([34|) . we then deduce that at C = the 
eigenvalues of onsite discrete dark solitons are given by 
= 4(A + 7)7, which becomes the continuous band 



for nonzero C, and Q 



(0) 



■7 



The continuation of the eigenvalue (^^^ for nonzero C 

can be calculated from Eq. (|35]) . The coefficient of rju"^ 
in this case is given by 



/:+(o)/:_(o)-i](°) = I 



0, 



4A7-A^ + 57^ riy^O. 



n = 0, . 
2 (40) 



The solvability condition for ([35]) then requires that /o = 
(4A - n^^^)vo^ = 0. Setting r]^^^ ^ 0, we deduce that 
Q'^^^ = —4 A. Hence the eigenvalue of an onsite dark 
soliton for small C is given by 



(41) 



Initially, i.e. for C = 0, the relative positions of the 
eigenvalue and the continuous spectrum can be divided 
into two cases, according to whether 7 ^ 7th = A/5. 
When C = and 7 < 7th (7 > 7th) the eigenvalue 
dHj) will be above (below) the continuous spectrum, as 
sketched in Fig. [3 These relative positions determine 
the instability mechanism for an onsite dark soliton, as 
we now describe. 

y<yth 



y>yth 



Dp 



1. Onsite dark solitons 



With errors of order , an onsite dark soliton is given 



by 



Un= \ 



( -y /AT^ , n = -2,-3,..., 

-V/AT7+ |C/v/AT7, n = -l, 

0,^ n = 0, 

V /AT7 - ^^/V'AT7, n = l, 

VA + 7, n = 2,3, .... 



For this configuration, 
£+(0)£_(0) = 



A2-72, n = 0, 
4(A + 7)7, HT^O. 



(38) 
(39) 



FIG. 7: As Fig. [T] but for a stable onsite dark soliton. 

For 7 < A/5, the instability is due to a colhsion be- 
tween the eigenvalue (|^T]) and which approximately 
occurs when 7 = 7^^., where 



-7A - ?C + i\/9A2 - 28AC - 16C2; 



(42) 



note that this critical value is meaningful only when C < 
9A/(14 + 2^85). For 7 > A/5, the instability is caused 
by the eigenvalue (|4T]) becoming negative, which occurs 
when 7 = 7^^, where 



4AC; 



(43) 
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this value is meaningful only when C < A/4. 

Furthermore, if we include terms up to O(C^), we ob- 
tain 

1] = - 7^ - 4AC + + 0{C^) (44) 

as the eigenvalue of an onsite discrete dark soliton. Using 
this expression, we find the critical value of 7 indicating 
the onset of instability to be 

7cr = -7A - ?G + iV9A2-28AC + 4C2, (45) 
5 5 5 

for 7 < 0.2A and 

-fl = VA2-4AC + 4C2, (46) 

for 7 > 0.2A. 



Hence the eigenvalues bifurcating from the lower edge of 
the continuous spectrum are given by 

= 4(A + 7)7 - (127 + 2A)C + (9(C'), (51) 

and 

= 4(A + 7)7 + (2A - 47)C + 0{C^). (52) 

y<0.5A 



2. Intersite modes 

Intersite discrete dark solitons are given, with errors of 
0(C2), by 



- VA+^ , n = -2,-3,..., 

-VT+T + C'/V^+T, n = -l, , 

n = 0, ^^'> 
VA + 7, n = l,2,.... 



Starting from Eq. (|34|) . we then find 
£+(0)£_(0)=4(A + 7)7 



(48) 



for all n, from which we deduce that there is a single 
leading-order eigenvalue, given by = 4(A + 7)7, 

with infinite multiplicity. This eigenvalue then expands 
to form the continuous spectrum for nonzero C. 

Because a localized structure must have an eigenvalue, 
we infer that an eigenvalue will bifurcate from the lower 
edge of the continuous spectrum. This bifurcating eigen- 
value may be calculated from Eq. (|35|) . Because 



£f (o)/:L'^(o)-i](o) = 



(49) 



for all n, the solvability condition for Eq. ([55|) requires 
/„ = for all n. A simple calculation then yields 

f =/ [4A + I67 + (2A + 47)A2 + r?i°\ n=-l,0, 
\ [(2A + 47)A2 + l^W]7?i°\ UT^-l.O. 



Taking r], 
equations 



(0) 



(50) 

for n 7^ —1,0 leaves the two nontrivial 



(87 + l^(i))7?i°] + (2A + 47)7?(°^ = 0, 
(87 + 1^(1) )7?r + (2A + A^)r,^°l = 0, 

from which we see thatr^i^^ =±Vo^' Thus we obtain two 
possibilities for the 0{C) contribution to the eigenvalue, 
given by 

n^^^ = -(127 + 2A), n^^^ = 2A - 47. 



y>03A 



FIG. 8: As Fig. [T] but for a stable intersite dark soliton. 

A simple analysis shows that 0.2 only when 7 > 

A/2. The sketch in Fig. [8] then illustrates that instability 
is caused by l^i becoming negative. This consideration 
gives the critical 7 as a function of the coupling constant 
C to be 



7c: 



-^A + + i - 4AC + 9C2. 



(53) 



When there are two eigenvalues (^^i and (^2), ^2 de- 
creases more slowly than l^i, in such a way that for 
7 > A/2 the instability is still caused by Vti becoming 
negative. 



B. Comparison with numerical computations 

1. Onsite dark solitons 

We now compare our analytical results with corre- 
sponding numerical calculations. As for bright solitons, 
for illustrative purposes we set A = 1. 

We start by checking the validity of our analytical ap- 
proximation for the critical eigenvalue associated with 
the phase mode. As explained above, the change in 
the position of the eigenvalues relative to the continu- 
ous spectrum at C = occurs at 7 = 1/5. Therefore we 
consider the two values 7 = 0.1 and 7 = 0.6, represent- 
ing both cases. Figure [9] depicts a comparison between 
our analytical result Eq. (|4T]) and the numerical compu- 
tations, from which we conclude that the prediction is 
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FIG. 9: Comparisons between the critical eigenvalue for on- 
site dark solitons obtained numerically (solid lines) and an- 
alytically using Eq. (|4T]) (dashed lines) for 7 = 0.1 (upper 
curves) and 7 = 0.6 (lower curves). An approximation that 
explicitly includes the next term in expansion Eq. (|44|) is also 
plotted (dotted lines). 



quite accurate for small C. The accuracy can be im- 
proved if one includes further orders in the perturbative 
expansion Eq. (|44|) . and this improvement is shown in 
the same figure by the dotted line. 

The eigenvalue structure of onsite dark solitons is de- 
picted in Fig. [TOl left and right panels refer respectively 
to conditions just before and just after a collision result- 
ing in an instability. 

As sketched in Fig. [71 for 7 < 1/5 the instability is 
caused by a collision between the eigenvalue and one edge 
of the continuous spectrum. On the other hand, when 
7 > 1/5 the instability is caused by a collision between 
the eigenvalue and its twin at the origin (see the bottom 
panels of Fig. [TQ|). 

We now proceed to evaluate the (in) stability region of 
this solution in (C, 7) space. Shown in Fig. [TT] is again 
the maximum of the imaginary part of the eigenvalue, to- 
gether with our approximation to the (in) stability bound- 
ary. The white dashed line represents Eq. (|43|) , corre- 
sponding to the instability caused by the collision with 
the continuous spectrum. Equation (|42]) is represented 
by the white dash-dotted line, which corresponds to the 
other instability mechanism. In addition, pink dashed 
and dash-dotted lines show, respectively, Eq. (|^ and 
Eq. (|45|) , where a better analytical approximation is ob- 
tained. 

An important observation from the figure is that there 
is an interval of values of 7 in which the onsite dark 
soliton is always stable, for any value of the coupling 
constant C. This indicates that a parametric driving can 
fully suppress the oscillatory instability reported for the 
first time in [21j . 



2. Intersite dark solitons 

Now we examine intersite dark solitons. 

Firstly, Fig. [12] shows the analytical prediction for the 
critical eigenvalue, given by Eq. ([5T]) , compared to nu- 
merical results. We see that the approximation is excel- 
lent for small C and its range of validity is wider for larger 
values of 7. The eigenvalue structure of this configuration 
is shown in Fig. [13] for the two values 7 = 0.1,0.8. The 
mechanism of instability explained in the section above 
can be seen clearly in the top panels of Fig. [13] 

It is interesting to note that a parametric driving 
can also fully suppress the oscillatory instability of an 
intersite dark soliton. As shown in the bottom pan- 
els of Fig. [131 there are values of the parameter 7 for 
which no instability-inducing collision ever occurs. The 
(in) stability region of this configuration is summarized 
in Fig. [HI where we see that for any C and 7 > 0.3 an 
intersite dark soliton is always stable. Our analytical pre- 
diction for the onset of instability is given by the dashed 
line in that figure. We observe that for relatively small 
C, the prediction of Eq. ([5T]) is reasonably close to the 
numerical results. 



IV. DISCUSSION 

In the sections above we discussed the existence and 
the stability of localized modes through our reduced 
DNLS equation ([5]). In this section, we confirm the 
relevance of our findings through solving numerically 
the original time-dependent equation ([T]). We use a 
Runge-Kutta integration method, with the initial con- 
dition (fri = 2eUn and (pn = 0, where Un is the static 
solution of the DNLS ([5|) and e is the small parameter 
of Sec. I. Throughout this section, we use the illustrative 
value e = 0.2. 

Shown in the left and right panels of Fig. [15] are the 
numerical evolution of a stable and unstable onsite bright 
soliton, respectively. From the right panel of the figure, 
we note that a parametric driving seems to destroy an 
unstable soliton. This observation is similar to the cor- 
responding observation for the dynamics of an unstable 
soliton in the DNLS equation ([5]) reported in Ref. [25|. 

In Fig. [16] we present the numerical evolution of in- 
tersite bright solitons for the same parameter values as 
those in Fig. [5l corresponding to each of the instability 
scenarios. From the panels in this figure, we see that 
the typical dynamics of the instability is in the form of 
soliton destruction or discharge of a traveling breather. 

We have also examined the dynamics of onsite dark 
solitons in the Klein-Gordon system ([1]). Shown in 
Fig. [T7] is the numerical evolution of a solution with the 
eigenvalue structure illustrated in Fig. [10] The instabil- 
ity of an unstable onsite dark soliton typically manifests 
itself in the form of oscillations in the location of the 
soliton center about its initial position (top right panel) 
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FIG. 10: The eigenvalue structure of on-site dark solitons for several values of 7 and C, as indicated in the caption of each 
panel. 
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FIG. 11: (Colour online) The (in)stability region of onsite 
dark solitons in the two-parameter (C, 7) space. The white 
and pink dashed lines respectively give the analytical approx- 
imations Eq. (|43)) and (|46)) . White and pink dash-dotted lines 
show Eqs. (|42|) and ([45^ : note that these curves are indistin- 
guishable in this plot. 
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FIG. 12: Comparisons between the critical eigenvalue for in- 
tersite dark solitons obtained numerically (solid lines) and 
analytically (dashed lines) using Eq. (|5T]) . for two values of 7. 
The upper curves correspond to 7 = 0.8, and the lower ones 
to 7 = 0.1. 



11 




3^ 
E 



0.1, C : 




0.1, C = 0.5 





(c)7 = 0.8, C = 0.5 (d)7 = 0.8, C = 2 

FIG. 13: The eigenvalue structure of intersite dark solitons with parameter values as indicated in the caption for each panel. 




FIG. 14: (Colour online) As Fig. [TT] but for an intersite dark 
soliton. The white dashed line is our analytical approximation 
Eq. (pl) . 



or oscillations in the width of the soliton (bottom right 
panel). 

Finally, we illustrate the dynamical behavior of an un- 



stable intersite dark soliton in Fig. [181 from which we 
see that the instability makes the soliton travel. This 
dynamics is similar to that reported in Ref. i20|. 



V. CONCLUSION 

In this paper, we have considered a parametrically 
driven Klein-Gordon system describing nanoelectrome- 
chanical systems. Using a multiscale expansion method 
we have reduced the system to a parametrically driven 
discrete nonlinear Schrodinger equation. Analytical and 
numerical calculations have been performed to deter- 
mine the existence and stability of fundamental bright 
and dark discrete solitons in the Klein-Gordon system 
through use of the Schrodinger equation. We have shown 
that the presence of a parametric driving can destabilize 
an onsite bright soliton. On the other hand, a parametric 
driving has also been shown to stabilize intersite bright 
and dark discrete solitons. We even found an interval 
in 7 for which a discrete dark soliton is stable for any 
value of the coupling constant, i.e. a parametric driving 
can suppress oscillatory instabilities. Stability windows 
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FIG. 15: (Colour online) The spatio-temporal evolution of an onsite bright soliton governed by the original time-dependent 
parametrically driven Klein-Gordon system ([T]), with e = 0.2 and 7 = 0.1. The left and right panels show a stable and unstable 
soliton, at C = 0.1 and C = 1, respectively. 



for all the fundamental solitons have been presented and 
approximations using perturbation theory have been de- 
rived to accompany the numerical results. Numerical 
integrations of the original Klein-Gordon system have 



demonstrated that our analytical and numerical inves- 
tigations of the discrete nonlinear Schrodinger equation 
provide a useful guide to behavior in the original system. 
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(a)7 = 0.1, C = 0.05 (b)7 = 0.1, C = 0.3 




FIG. 16: (Colour online) As Fig. I15[ but for an intersite bright soliton, with parameter values as indicated in the caption for 
each panel. The initial profile in each panel corresponds to the same parameters as in Fig. [S] 
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(c)7 = 0.6, C = 0.01 (d)7 = 0.6, C = 1 



FIG. 17: (Colour online) As Fig. [151 but for on-site dark solitons. The parameter values are as in Fig. [TOl 




FIG. 18: (Colour online) As Fig. 1151 but for an intersite dark soliton with 7 = 0.1. The left panel shows the evolution of a 
stable dark soliton with C — 0.05, while the right panel shows the evolution of an unstable dark soliton with C — 0.5. 



